home *** CD-ROM | disk | FTP | other *** search
/ Linux Cubed Series 2: Applications / Linux Cubed Series 2 - Applications.iso / math / gle-3.000 / gle-3 / gle / lib / feyn.gle < prev   
Text File  |  1995-03-08  |  51KB  |  1,627 lines

  1. ! feyn.gle - drawing Feynman diagrams, version 1.0 
  2. ! Copyright (c) Andrey Grozin, 1993-1995 
  3. ! variable names beginning with F are reserved for internal use 
  4.  
  5. !         straight line | arc 
  6. ! ---------------------------- 
  7. ! fr    - -1            | radius 
  8. ! fx,fy - initial point | center 
  9. ! fa    - slope         | slope: center->initial 
  10. ! fl    - length        | angle length 
  11. ! fc    -               | +1 - counterclockwise 
  12. ! ffx,ffy         - final point 
  13. ! fx1,fy1,...     - temporary coordinates 
  14. ! f1,f2,...       - temporary variables 
  15. ! fi,fm,fd        - for loops 
  16.  
  17.  
  18. pi = 3.14159265358979323846264338328 
  19. ! sin table 
  20. fs1 = 0.156434465040230869010105319467 
  21. fs2 = 0.309016994374947424102293417183 
  22. fs3 = 0.453990499739546791560408366358 
  23. fs4 = 0.587785252292473129168705954639 
  24. fs5 = 0.707106781186547524400844362105 
  25. fs6 = 0.809016994374947424102293417183 
  26. fs7 = 0.891006524188367862359709571414 
  27. fs8 = 0.951056516295153572116439333379 
  28. fs9 = 0.987688340595137726190040247693 
  29.  
  30. fe = 0.01 ! lines shorter than FE cm are not drawn 
  31.  
  32. set lwidth 0.05 ! line width 
  33. set cap  round  ! end of line 
  34. set join round  ! join lines 
  35. set hei 0.5     ! chaaracter height 
  36. set font texmi  ! PostScript italic 
  37. set just cc     ! centered 
  38.  
  39. ! simple debugging tool 
  40. DebugX = 1 
  41. DebugY = 1 
  42. DebugD = 1 
  43. Debug = 0 
  44. sub fdebug l t$ v 
  45.     if Debug>l-0.5 then 
  46.         gsave 
  47.         amove DebugX DebugY 
  48.         write t$+num$(v) 
  49.         grestore 
  50.     DebugY = DebugY+DebugD 
  51.     end if 
  52. end sub 
  53.  
  54. ! internal - setup for straight line 
  55. sub feyn x y 
  56.     ffx = x 
  57.     ffy = y 
  58.     fr = -1 
  59.     fx = xpos() 
  60.     fy = ypos() 
  61.     x = x-fx 
  62.     y = y-fy 
  63.     fl = sqrt(sqr(x)+sqr(y)) 
  64.     if fl<fe then 
  65.         fl = -1 
  66.     else ! finding slope 
  67.         if abs(y)<abs(x) then 
  68.             if x>0 then 
  69.                 fa = 180*atn(y/x)/pi 
  70.             else 
  71.                 fa = 180*(1+atn(y/x)/pi) 
  72.             end if 
  73.         else 
  74.             if y>0 then 
  75.                 fa = 180*(0.5-atn(x/y)/pi) 
  76.             else 
  77.                 fa = 180*(1.5-atn(x/y)/pi) 
  78.             end if 
  79.         end if 
  80.     end if 
  81. end sub 
  82.  
  83. ! internal - setup for arc 
  84. sub feyn2 xx yy x y 
  85.     ffx = x 
  86.     ffy = y 
  87.     fx1 = xpos()-xx 
  88.     fy1 = ypos()-yy 
  89.     fx2 = x-xx 
  90.     fy2 = y-yy 
  91.     fr = fx1*fy2-fx2*fy1 
  92.     if abs(fr)<sqr(fe) then 
  93.         fr = -1 
  94.     else ! finding center 
  95.         fc = -sgn(fr) 
  96.         fx = xx+0.5/fr*(sqr(fy1)*fy2-sqr(fy2)*fy1+sqr(fx1)*fy2-sqr(fx2)*fy1) 
  97.         fy = yy+0.5/fr*(sqr(fx2)*fx1-sqr(fx1)*fx2+sqr(fy2)*fx1-sqr(fy1)*fx2) 
  98.         fx1 = xpos()-fx 
  99.         fy1 = ypos()-fy 
  100.         fr = sqrt(sqr(fx1)+sqr(fy1)) 
  101.         ! finding slope center - initial 
  102.         if abs(fy1)<abs(fx1) then 
  103.             if fx1>0 then 
  104.                 fa = 180*atn(fy1/fx1)/pi 
  105.             else 
  106.                 fa = 180*(1+atn(fy1/fx1)/pi) 
  107.             end if 
  108.         else 
  109.             if fy1>0 then 
  110.                 fa = 180*(0.5-atn(fx1/fy1)/pi) 
  111.             else 
  112.                 fa = 180*(1.5-atn(fx1/fy1)/pi) 
  113.             end if 
  114.         end if 
  115.         ! finding angular length 
  116.         fx2 = x-fx 
  117.         fy2 = y-fy 
  118.         f1 = (fx1*fy2-fx2*fy1)/sqr(fr)*fc ! sin 
  119.         f2 = (fx1*fx2+fy1*fy2)/sqr(fr)    ! cos 
  120.         if f1<0 then                      ! pi .. 2*pi 
  121.             if abs(f1)>abs(f2) then       ! 5*pi/4 .. 7*pi/4 
  122.                 if f2<0 then              ! 5*pi/4 .. 3*pi/2 
  123.                     fl = 180*(1.5-atn(f2/f1)/pi) 
  124.                 else                      ! 3*pi/2 .. 7*pi/4 
  125.                     fl = 180*(1.5+atn(-f2/f1)/pi) 
  126.                 end if 
  127.             else                          ! pi .. 5*pi/4 or 7*pi/4 .. 2*pi 
  128.                 if f2<0 then              ! pi .. 5*pi/4 
  129.                     fl = 180*(1+atn(f1/f2)/pi) 
  130.                 else                      ! 7*pi/4 .. 2*pi 
  131.                     fl = 180*(2-atn(-f1/f2)/pi) 
  132.                 end if 
  133.             end if 
  134.         else                              ! 0 .. pi 
  135.             if abs(f1)>abs(f2) then       ! pi/4 .. 3*pi/4 
  136.                 if f2<0 then              ! pi/2 .. 3*pi/4 
  137.                     fl = 180*(0.5+atn(-f2/f1)/pi) 
  138.                 else                      ! pi/4 .. pi/2 
  139.                     fl = 180*(0.5-atn(f2/f1)/pi) 
  140.                 end if 
  141.             else                          ! 0 .. pi/4 or 3*pi/4 .. pi 
  142.                 if f2<0 then              ! 3*pi/4 .. pi 
  143.                     fl = 180*(1-atn(-f1/f2)/pi) 
  144.                 else                      ! 0 .. pi/4 
  145.                     fl = 180*atn(f1/f2)/pi 
  146.                 end if 
  147.             end if 
  148.         end if 
  149.     end if 
  150. end sub 
  151.  
  152. ! Fermion line from the current point to x,y 
  153. sub Fermion x y 
  154.     @feyn x y 
  155.     if fl>0 then 
  156.         begin origin 
  157.             begin rotate fa 
  158.                 aline fl 0 
  159.             end rotate 
  160.         end origin 
  161.     amove ffx ffy 
  162.     end if 
  163. end sub 
  164.  
  165. ! Fermion line from the current point via xx,yy to x,y 
  166. sub Fermion2 xx yy x y 
  167.     @feyn2 xx yy x y 
  168.     if fr<0 then 
  169.         @Fermion x y 
  170.     else 
  171.         amove fx fy 
  172.         begin origin 
  173.             if fc>0 then 
  174.                 arc fr fa fa+fl 
  175.             else 
  176.                 arc fr fa-fl fa 
  177.             end if 
  178.         end origin 
  179.     amove ffx ffy 
  180.     end if 
  181. end sub 
  182.  
  183. DoubleA = 0.05 ! half distance between lines 
  184.  
  185. ! Double line from the current point to x,y 
  186. sub Double x y 
  187.     @feyn x y 
  188.     if fl>0 then 
  189.         begin origin 
  190.             begin rotate fa 
  191.         amove 0   DoubleA 
  192.         aline fl  DoubleA 
  193.         amove 0  -DoubleA 
  194.         aline fl -DoubleA 
  195.             end rotate 
  196.         end origin 
  197.     amove ffx ffy 
  198.     end if 
  199. end sub 
  200.  
  201. ! Double line from the current point via xx,yy to x,y 
  202. sub Double2 xx yy x y 
  203.     @feyn2 xx yy x y 
  204.     if fr<DoubleA+fe then 
  205.         @Double x y 
  206.     else 
  207.         amove fx fy 
  208.         begin origin 
  209.             if fc>0 then 
  210.                 arc fr+DoubleA fa fa+fl 
  211.         arc fr-DoubleA fa fa+fl 
  212.             else 
  213.                 arc fr+DoubleA fa-fl fa 
  214.         arc fr-DoubleA fa-fl fa 
  215.             end if 
  216.         end origin 
  217.     amove ffx ffy 
  218.     end if 
  219. end sub 
  220.  
  221. DashL = 0.3  ! period 
  222. DashF = 0.5  ! fraction of the period that is filled 
  223.  
  224. ! Higgs (dashed) line from the current point to x,y 
  225. sub Dash x y 
  226.     @feyn x y 
  227.     if fl>0 then 
  228.         fn = fix(fl/DashL+1.5-DashF) 
  229.         if fn<1 then 
  230.             fn = 1 
  231.         end if 
  232.         fd = fl/(fn-1+DashF) 
  233.         fm = fl-0.5*DashF*fd 
  234.         f1 = fd*DashF 
  235.         begin origin 
  236.             begin rotate fa 
  237.                 for fi = 0 to fm step fd 
  238.                     begin translate fi 0 
  239.                         amove 0 0 
  240.                         aline f1 0 
  241.                     end translate 
  242.                 next fi 
  243.             end rotate 
  244.         end origin 
  245.     amove ffx ffy 
  246.     end if 
  247. end sub 
  248.  
  249. ! Higgs (dashed) line from the current point via xx,yy to x,y 
  250. sub Dash2 xx yy x y 
  251.     @feyn2 xx yy x y 
  252.     if fr<0 then 
  253.         @Dash x y 
  254.     else 
  255.         fn = fix(fr*fl*pi/180/DashL+1.5-DashF) 
  256.         if fn<1 then 
  257.             fn = 1 
  258.         end if 
  259.         fd = fl/(fn-1+DashF) 
  260.         fm = fl-0.5*DashF*fd 
  261.         f1 = fd*DashF 
  262.         amove fx fy 
  263.         begin origin 
  264.             if fc>0 then 
  265.                 for fi = 0 to fm step fd 
  266.                     arc fr fa+fi fa+fi+f1 
  267.                 next fi 
  268.             else 
  269.                 for fi = 0 to fm step fd 
  270.                     narc fr fa-fi fa-fi-f1 
  271.                 next fi 
  272.             end if 
  273.         end origin 
  274.     amove ffx ffy 
  275.     end if 
  276. end sub 
  277.  
  278. PhotonL = 0.1  ! half-wavelength 
  279. PhotonA = 0.1  ! amplitude 
  280. PhotonN = 0    ! correction to the number of half-waves 
  281.  
  282. ! Vector boson (zigzag) line from the current point to x,y 
  283. sub Zigzag x y 
  284.     @feyn x y 
  285.     if fl>0 then 
  286.         fn = fix(fl/PhotonL+0.5+PhotonN) 
  287.         if fn<1 then 
  288.             fn = 1 
  289.         end if 
  290.     PhotonN = 0 
  291.         fd = fl/fn 
  292.         fm = fl-0.5*fd 
  293.         f1 = PhotonA 
  294.         begin origin 
  295.             begin rotate fa 
  296.                 for fi = 0 to fm step fd 
  297.                     begin translate fi 0 
  298.             amove 0 0 
  299.             aline 0.5*fd f1 
  300.             aline fd 0 
  301.                     end translate 
  302.                     f1 = -f1 
  303.                 next fi 
  304.             end rotate 
  305.         end origin 
  306.     amove ffx ffy 
  307.     end if 
  308. end sub 
  309.  
  310. ! Vector boson (zigzag) line from the current point via xx,yy to x,y 
  311. sub Zigzag2 xx yy x y 
  312.     @feyn2 xx yy x y 
  313.     if fr<PhotonA+fe then 
  314.         @Zigzag x y 
  315.     else 
  316.         fn = fix(fr*fl*pi/180/PhotonL+0.5+PhotonN) 
  317.         if fn<1 then 
  318.             fn = 1 
  319.         end if 
  320.     PhotonN = 0 
  321.         fd = fl/fn 
  322.     f1 = fc*PhotonA/fr 
  323.     f2 = sin(pi*fd/360) 
  324.     f3 = cos(pi*fd/360) 
  325.     f4 = sqrt(f2^4+f1^2*f3^2*(1+2*f2^2)-f1^4*f2^2*f3^2) 
  326.     if f1>0 then 
  327.         f2 = (1+f1)*f2*f3*(1-(f1-f4)/(f2^2+f1^2*f3^2)) 
  328.     else 
  329.         f2 = (1+f1)*f2*f3*(1-(f1+f4)/(f2^2+f1^2*f3^2)) 
  330.     end if 
  331.     f2 = 180/pi*atn(f2/sqrt(1-f2^2)) 
  332.     if abs(fn-2*fix(0.5*fn))>0.1 then 
  333.         if abs(fn-1)<0.1 then 
  334.             fd = fl 
  335.         f2 = 0.5*fl 
  336.         else 
  337.             for fi=1 to 10 
  338.             fd = (fl-2*f2)/(fn-1) 
  339.             f2 = sin(pi*fd/360) 
  340.             f3 = cos(pi*fd/360) 
  341.             f4 = sqrt(f2^4+f1^2*f3^2*(1+2*f2^2)-f1^4*f2^2*f3^2) 
  342.             if f1>0 then 
  343.                 f2 = (1+f1)*f2*f3*(1-(f1-f4)/(f2^2+f1^2*f3^2)) 
  344.             else 
  345.                 f2 = (1+f1)*f2*f3*(1-(f1+f4)/(f2^2+f1^2*f3^2)) 
  346.             end if 
  347.             f2 = 180/pi*atn(f2/sqrt(1-f2^2)) 
  348.         next fi 
  349.         end if 
  350.     end if 
  351.     fm = f2+(2*fix(0.5*(fn+1)+0.1)-3)*fd 
  352.     fx1 = fr*(1+f1)*cos(pi*fd/180) 
  353.     fy1 = fr*(1+f1)*sin(pi*fd/180) 
  354.     fx2 = fr*(1-f1)*cos(pi*fd/90) 
  355.     fy2 = fr*(1-f1)*sin(pi*fd/90) 
  356.     fx3 = fr*(1-f1)*cos(pi*f2/180) 
  357.     fy3 = fr*(1-f1)*sin(pi*f2/180) 
  358.         amove fx fy 
  359.         begin origin 
  360.             if fc>0 then 
  361.         begin rotate fa 
  362.             amove fr 0 
  363.             aline fx3 fy3 
  364.         end rotate 
  365.         if fn>2.5 then 
  366.             for fi = f2 to fm step 2*fd 
  367.                 begin rotate fa+fi 
  368.                 amove fr*(1-f1) 0 
  369.                 aline fx1 fy1 
  370.                 aline fx2 fy2 
  371.             end rotate 
  372.             next fi 
  373.         end if 
  374.         f2 = fl-fm-fd 
  375.         if abs(fn-2*fix(0.5*fn))<0.1 then 
  376.             begin rotate fa+fm+fd 
  377.             amove fr*(1-f1) 0 
  378.             aline fx1 fy1 
  379.             end rotate 
  380.             f2 = f2-fd 
  381.             f1 = -f1 
  382.         end if 
  383.         fx3 = fr*(1-f1)*cos(pi*f2/180) 
  384.         fy3 = fr*(1-f1)*sin(pi*f2/180) 
  385.         begin rotate fa+fl 
  386.             amove fx3 -fy3 
  387.             aline fr 0 
  388.         end rotate 
  389.             else 
  390.         begin rotate fa 
  391.             amove fr 0 
  392.             aline fx3 -fy3 
  393.         end rotate 
  394.         if fn>2.5 then 
  395.             for fi = f2 to fm step 2*fd 
  396.                 begin rotate fa-fi 
  397.                 amove fr*(1-f1) 0 
  398.                 aline fx1 -fy1 
  399.                 aline fx2 -fy2 
  400.             end rotate 
  401.             next fi 
  402.         end if 
  403.         f2 = fl-fm-fd 
  404.         if abs(fn-2*fix(0.5*fn))<0.1 then 
  405.             begin rotate fa-fm-fd 
  406.             amove fr*(1-f1) 0 
  407.             aline fx1 -fy1 
  408.             end rotate 
  409.             f2 = f2-fd 
  410.             f1 = -f1 
  411.         end if 
  412.         fx3 = fr*(1-f1)*cos(pi*f2/180) 
  413.         fy3 = fr*(1-f1)*sin(pi*f2/180) 
  414.         begin rotate fa-fl 
  415.             amove fx3 fy3 
  416.             aline fr 0 
  417.         end rotate 
  418.             end if 
  419.         end origin 
  420.     amove ffx ffy 
  421.     end if 
  422. end sub 
  423.  
  424. ! Photon line from the current point to x,y 
  425. sub Photon x y 
  426.     @feyn x y 
  427.     if fl>0 then 
  428.         fn = fix(fl/PhotonL+0.5+PhotonN) 
  429.         if fn<1 then 
  430.             fn = 1 
  431.         end if 
  432.     PhotonN = 0 
  433.         fd = fl/fn 
  434.         fm = fl-0.5*fd 
  435.         f1 = PhotonA 
  436.         begin origin 
  437.             begin rotate fa 
  438.                 for fi = 0 to fm step fd 
  439.                     begin translate fi 0 
  440.             amove 0       0 
  441.             aline 0.05*fd fs1*f1 
  442.             aline 0.1*fd  fs2*f1 
  443.             aline 0.15*fd fs3*f1 
  444.             aline 0.2*fd  fs4*f1 
  445.             aline 0.25*fd fs5*f1 
  446.             aline 0.3*fd  fs6*f1 
  447.             aline 0.35*fd fs7*f1 
  448.             aline 0.4*fd  fs8*f1 
  449.             aline 0.45*fd fs9*f1 
  450.             aline 0.5*fd  f1 
  451.             aline 0.55*fd fs9*f1 
  452.             aline 0.6*fd  fs8*f1 
  453.             aline 0.65*fd fs7*f1 
  454.             aline 0.7*fd  fs6*f1 
  455.             aline 0.75*fd fs5*f1 
  456.             aline 0.8*fd  fs4*f1 
  457.             aline 0.85*fd fs3*f1 
  458.             aline 0.9*fd  fs2*f1 
  459.             aline 0.95*fd fs1*f1 
  460.             aline fd      0 
  461.                     end translate 
  462.                     f1 = -f1 
  463.                 next fi 
  464.             end rotate 
  465.         end origin 
  466.     amove ffx ffy 
  467.     end if 
  468. end sub 
  469.  
  470. ! Photon line from the current point via xx,yy to x,y 
  471. sub Photon2 xx yy x y 
  472.     @feyn2 xx yy x y 
  473.     if fr<PhotonA+fe then 
  474.         @Photon x y 
  475.     else 
  476.         fn = fix(fr*fl*pi/180/PhotonL+0.5+PhotonN) 
  477.         if fn<1 then 
  478.             fn = 1 
  479.         end if 
  480.     PhotonN = 0 
  481.         fd = fl/fn 
  482.     f1 = fc*PhotonA/fr 
  483.     f2 = sin(pi*fd/360) 
  484.     f3 = cos(pi*fd/360) 
  485.     f4 = sqrt(f2^4+f1^2*f3^2*(1+2*f2^2)-f1^4*f2^2*f3^2) 
  486.     if f1>0 then 
  487.         f2 = (1+f1)*f2*f3*(1-(f1-f4)/(f2^2+f1^2*f3^2)) 
  488.     else 
  489.         f2 = (1+f1)*f2*f3*(1-(f1+f4)/(f2^2+f1^2*f3^2)) 
  490.     end if 
  491.     f2 = 180/pi*atn(f2/sqrt(1-f2^2)) 
  492.     if abs(fn-2*fix(0.5*fn))>0.1 then 
  493.         if abs(fn-1)<0.1 then 
  494.             fd = fl 
  495.         f2 = 0.5*fl 
  496.         else 
  497.             for fi=1 to 10 
  498.             fd = (fl-2*f2)/(fn-1) 
  499.             f2 = sin(pi*fd/360) 
  500.             f3 = cos(pi*fd/360) 
  501.             f4 = sqrt(f2^4+f1^2*f3^2*(1+2*f2^2)-f1^4*f2^2*f3^2) 
  502.             if f1>0 then 
  503.                 f2 = (1+f1)*f2*f3*(1-(f1-f4)/(f2^2+f1^2*f3^2)) 
  504.             else 
  505.                 f2 = (1+f1)*f2*f3*(1-(f1+f4)/(f2^2+f1^2*f3^2)) 
  506.             end if 
  507.             f2 = 180/pi*atn(f2/sqrt(1-f2^2)) 
  508.         next fi 
  509.         end if 
  510.     end if 
  511.     fm = f2+(2*fix(0.5*(fn+1)+0.1)-3)*fd 
  512.         amove fx fy 
  513.         begin origin 
  514.             if fc>0 then 
  515.         begin rotate fa 
  516.             amove fr 0 
  517.             fx1 = fr*(1-f1*fs1)*cos(1*pi*f2/1800) 
  518.             fy1 = fr*(1-f1*fs1)*sin(1*pi*f2/1800) 
  519.             aline fx1 fy1 
  520.                 fx1 = fr*(1-f1*fs2)*cos(2*pi*f2/1800) 
  521.             fy1 = fr*(1-f1*fs2)*sin(2*pi*f2/1800) 
  522.             aline fx1 fy1 
  523.                 fx1 = fr*(1-f1*fs3)*cos(3*pi*f2/1800) 
  524.             fy1 = fr*(1-f1*fs3)*sin(3*pi*f2/1800) 
  525.             aline fx1 fy1 
  526.                 fx1 = fr*(1-f1*fs4)*cos(4*pi*f2/1800) 
  527.             fy1 = fr*(1-f1*fs4)*sin(4*pi*f2/1800) 
  528.             aline fx1 fy1 
  529.                 fx1 = fr*(1-f1*fs5)*cos(5*pi*f2/1800) 
  530.             fy1 = fr*(1-f1*fs5)*sin(5*pi*f2/1800) 
  531.             aline fx1 fy1 
  532.                 fx1 = fr*(1-f1*fs6)*cos(6*pi*f2/1800) 
  533.             fy1 = fr*(1-f1*fs6)*sin(6*pi*f2/1800) 
  534.             aline fx1 fy1 
  535.                 fx1 = fr*(1-f1*fs7)*cos(7*pi*f2/1800) 
  536.             fy1 = fr*(1-f1*fs7)*sin(7*pi*f2/1800) 
  537.             aline fx1 fy1 
  538.                 fx1 = fr*(1-f1*fs8)*cos(8*pi*f2/1800) 
  539.             fy1 = fr*(1-f1*fs8)*sin(8*pi*f2/1800) 
  540.             aline fx1 fy1 
  541.                 fx1 = fr*(1-f1*fs9)*cos(9*pi*f2/1800) 
  542.             fy1 = fr*(1-f1*fs9)*sin(9*pi*f2/1800) 
  543.             aline fx1 fy1 
  544.                 fx1 = fr*(1-f1)*cos(pi*f2/180) 
  545.             fy1 = fr*(1-f1)*sin(pi*f2/180) 
  546.             aline fx1 fy1 
  547.         end rotate 
  548.         if fn>2.5 then 
  549.             for fi = f2 to fm step 2*fd 
  550.                 begin rotate fa+fi 
  551.                 amove fr*(1-f1) 0 
  552.                 fx1 = fr*(1-f1*fs9)*cos(1*pi*fd/3600) 
  553.                 fy1 = fr*(1-f1*fs9)*sin(1*pi*fd/3600) 
  554.                 aline fx1 fy1 
  555.                 fx1 = fr*(1-f1*fs8)*cos(2*pi*fd/3600) 
  556.                 fy1 = fr*(1-f1*fs8)*sin(2*pi*fd/3600) 
  557.                 aline fx1 fy1 
  558.                 fx1 = fr*(1-f1*fs7)*cos(3*pi*fd/3600) 
  559.                 fy1 = fr*(1-f1*fs7)*sin(3*pi*fd/3600) 
  560.                 aline fx1 fy1 
  561.                 fx1 = fr*(1-f1*fs6)*cos(4*pi*fd/3600) 
  562.                 fy1 = fr*(1-f1*fs6)*sin(4*pi*fd/3600) 
  563.                 aline fx1 fy1 
  564.                 fx1 = fr*(1-f1*fs5)*cos(5*pi*fd/3600) 
  565.                 fy1 = fr*(1-f1*fs5)*sin(5*pi*fd/3600) 
  566.                 aline fx1 fy1 
  567.                 fx1 = fr*(1-f1*fs4)*cos(6*pi*fd/3600) 
  568.                 fy1 = fr*(1-f1*fs4)*sin(6*pi*fd/3600) 
  569.                 aline fx1 fy1 
  570.                 fx1 = fr*(1-f1*fs3)*cos(7*pi*fd/3600) 
  571.                 fy1 = fr*(1-f1*fs3)*sin(7*pi*fd/3600) 
  572.                 aline fx1 fy1 
  573.                 fx1 = fr*(1-f1*fs2)*cos(8*pi*fd/3600) 
  574.                 fy1 = fr*(1-f1*fs2)*sin(8*pi*fd/3600) 
  575.                 aline fx1 fy1 
  576.                 fx1 = fr*(1-f1*fs1)*cos(9*pi*fd/3600) 
  577.                 fy1 = fr*(1-f1*fs1)*sin(9*pi*fd/3600) 
  578.                 aline fx1 fy1 
  579.                 fx1 = fr*cos(pi*fd/360) 
  580.                 fy1 = fr*sin(pi*fd/360) 
  581.                 aline fx1 fy1 
  582.                 fx1 = fr*(1+f1*fs1)*cos(11*pi*fd/3600) 
  583.                 fy1 = fr*(1+f1*fs1)*sin(11*pi*fd/3600) 
  584.                 aline fx1 fy1 
  585.                 fx1 = fr*(1+f1*fs2)*cos(12*pi*fd/3600) 
  586.                 fy1 = fr*(1+f1*fs2)*sin(12*pi*fd/3600) 
  587.                 aline fx1 fy1 
  588.                 fx1 = fr*(1+f1*fs3)*cos(13*pi*fd/3600) 
  589.                 fy1 = fr*(1+f1*fs3)*sin(13*pi*fd/3600) 
  590.                 aline fx1 fy1 
  591.                 fx1 = fr*(1+f1*fs4)*cos(14*pi*fd/3600) 
  592.                 fy1 = fr*(1+f1*fs4)*sin(14*pi*fd/3600) 
  593.                 aline fx1 fy1 
  594.                 fx1 = fr*(1+f1*fs5)*cos(15*pi*fd/3600) 
  595.                 fy1 = fr*(1+f1*fs5)*sin(15*pi*fd/3600) 
  596.                 aline fx1 fy1 
  597.                 fx1 = fr*(1+f1*fs6)*cos(16*pi*fd/3600) 
  598.                 fy1 = fr*(1+f1*fs6)*sin(16*pi*fd/3600) 
  599.                 aline fx1 fy1 
  600.                 fx1 = fr*(1+f1*fs7)*cos(17*pi*fd/3600) 
  601.                 fy1 = fr*(1+f1*fs7)*sin(17*pi*fd/3600) 
  602.                 aline fx1 fy1 
  603.                 fx1 = fr*(1+f1*fs8)*cos(18*pi*fd/3600) 
  604.                 fy1 = fr*(1+f1*fs8)*sin(18*pi*fd/3600) 
  605.                 aline fx1 fy1 
  606.                 fx1 = fr*(1+f1*fs9)*cos(19*pi*fd/3600) 
  607.                 fy1 = fr*(1+f1*fs9)*sin(19*pi*fd/3600) 
  608.                 aline fx1 fy1 
  609.                 fx1 = fr*(1+f1)*cos(pi*fd/180) 
  610.                 fy1 = fr*(1+f1)*sin(pi*fd/180) 
  611.                 aline fx1 fy1 
  612.                 fx1 = fr*(1+f1*fs9)*cos(21*pi*fd/3600) 
  613.                 fy1 = fr*(1+f1*fs9)*sin(21*pi*fd/3600) 
  614.                 aline fx1 fy1 
  615.                 fx1 = fr*(1+f1*fs8)*cos(22*pi*fd/3600) 
  616.                 fy1 = fr*(1+f1*fs8)*sin(22*pi*fd/3600) 
  617.                 aline fx1 fy1 
  618.                 fx1 = fr*(1+f1*fs7)*cos(23*pi*fd/3600) 
  619.                 fy1 = fr*(1+f1*fs7)*sin(23*pi*fd/3600) 
  620.                 aline fx1 fy1 
  621.                 fx1 = fr*(1+f1*fs6)*cos(24*pi*fd/3600) 
  622.                 fy1 = fr*(1+f1*fs6)*sin(24*pi*fd/3600) 
  623.                 aline fx1 fy1 
  624.                 fx1 = fr*(1+f1*fs5)*cos(25*pi*fd/3600) 
  625.                 fy1 = fr*(1+f1*fs5)*sin(25*pi*fd/3600) 
  626.                 aline fx1 fy1 
  627.                 fx1 = fr*(1+f1*fs4)*cos(26*pi*fd/3600) 
  628.                 fy1 = fr*(1+f1*fs4)*sin(26*pi*fd/3600) 
  629.                 aline fx1 fy1 
  630.                 fx1 = fr*(1+f1*fs3)*cos(27*pi*fd/3600) 
  631.                 fy1 = fr*(1+f1*fs3)*sin(27*pi*fd/3600) 
  632.                 aline fx1 fy1 
  633.                 fx1 = fr*(1+f1*fs2)*cos(28*pi*fd/3600) 
  634.                 fy1 = fr*(1+f1*fs2)*sin(28*pi*fd/3600) 
  635.                 aline fx1 fy1 
  636.                 fx1 = fr*(1+f1*fs1)*cos(29*pi*fd/3600) 
  637.                 fy1 = fr*(1+f1*fs1)*sin(29*pi*fd/3600) 
  638.                 aline fx1 fy1 
  639.                 fx1 = fr*cos(pi*fd/120) 
  640.                 fy1 = fr*sin(pi*fd/120) 
  641.                 aline fx1 fy1 
  642.                 fx1 = fr*(1-f1*fs1)*cos(31*pi*fd/3600) 
  643.                 fy1 = fr*(1-f1*fs1)*sin(31*pi*fd/3600) 
  644.                 aline fx1 fy1 
  645.                 fx1 = fr*(1-f1*fs2)*cos(32*pi*fd/3600) 
  646.                 fy1 = fr*(1-f1*fs2)*sin(32*pi*fd/3600) 
  647.                 aline fx1 fy1 
  648.                 fx1 = fr*(1-f1*fs3)*cos(33*pi*fd/3600) 
  649.                 fy1 = fr*(1-f1*fs3)*sin(33*pi*fd/3600) 
  650.                 aline fx1 fy1 
  651.                 fx1 = fr*(1-f1*fs4)*cos(34*pi*fd/3600) 
  652.                 fy1 = fr*(1-f1*fs4)*sin(34*pi*fd/3600) 
  653.                 aline fx1 fy1 
  654.                 fx1 = fr*(1-f1*fs5)*cos(35*pi*fd/3600) 
  655.                 fy1 = fr*(1-f1*fs5)*sin(35*pi*fd/3600) 
  656.                 aline fx1 fy1 
  657.                 fx1 = fr*(1-f1*fs6)*cos(36*pi*fd/3600) 
  658.                 fy1 = fr*(1-f1*fs6)*sin(36*pi*fd/3600) 
  659.                 aline fx1 fy1 
  660.                 fx1 = fr*(1-f1*fs7)*cos(37*pi*fd/3600) 
  661.                 fy1 = fr*(1-f1*fs7)*sin(37*pi*fd/3600) 
  662.                 aline fx1 fy1 
  663.                 fx1 = fr*(1-f1*fs8)*cos(38*pi*fd/3600) 
  664.                 fy1 = fr*(1-f1*fs8)*sin(38*pi*fd/3600) 
  665.                 aline fx1 fy1 
  666.                 fx1 = fr*(1-f1*fs9)*cos(39*pi*fd/3600) 
  667.                 fy1 = fr*(1-f1*fs9)*sin(39*pi*fd/3600) 
  668.                 aline fx1 fy1 
  669.                 fx1 = fr*(1-f1)*cos(pi*fd/90) 
  670.                 fy1 = fr*(1-f1)*sin(pi*fd/90) 
  671.                 aline fx1 fy1 
  672.             end rotate 
  673.             next fi 
  674.         end if 
  675.         f2 = fl-fm-fd 
  676.         if abs(fn-2*fix(0.5*fn))<0.1 then 
  677.             begin rotate fa+fm+fd 
  678.             amove fr*(1-f1) 0 
  679.             fx1 = fr*(1-f1*fs9)*cos(1*pi*fd/3600) 
  680.             fy1 = fr*(1-f1*fs9)*sin(1*pi*fd/3600) 
  681.             aline fx1 fy1 
  682.             fx1 = fr*(1-f1*fs8)*cos(2*pi*fd/3600) 
  683.             fy1 = fr*(1-f1*fs8)*sin(2*pi*fd/3600) 
  684.             aline fx1 fy1 
  685.             fx1 = fr*(1-f1*fs7)*cos(3*pi*fd/3600) 
  686.             fy1 = fr*(1-f1*fs7)*sin(3*pi*fd/3600) 
  687.             aline fx1 fy1 
  688.             fx1 = fr*(1-f1*fs6)*cos(4*pi*fd/3600) 
  689.             fy1 = fr*(1-f1*fs6)*sin(4*pi*fd/3600) 
  690.             aline fx1 fy1 
  691.             fx1 = fr*(1-f1*fs5)*cos(5*pi*fd/3600) 
  692.             fy1 = fr*(1-f1*fs5)*sin(5*pi*fd/3600) 
  693.             aline fx1 fy1 
  694.             fx1 = fr*(1-f1*fs4)*cos(6*pi*fd/3600) 
  695.             fy1 = fr*(1-f1*fs4)*sin(6*pi*fd/3600) 
  696.             aline fx1 fy1 
  697.             fx1 = fr*(1-f1*fs3)*cos(7*pi*fd/3600) 
  698.             fy1 = fr*(1-f1*fs3)*sin(7*pi*fd/3600) 
  699.             aline fx1 fy1 
  700.             fx1 = fr*(1-f1*fs2)*cos(8*pi*fd/3600) 
  701.             fy1 = fr*(1-f1*fs2)*sin(8*pi*fd/3600) 
  702.             aline fx1 fy1 
  703.             fx1 = fr*(1-f1*fs1)*cos(9*pi*fd/3600) 
  704.             fy1 = fr*(1-f1*fs1)*sin(9*pi*fd/3600) 
  705.             aline fx1 fy1 
  706.             fx1 = fr*cos(pi*fd/360) 
  707.             fy1 = fr*sin(pi*fd/360) 
  708.             aline fx1 fy1 
  709.             fx1 = fr*(1+f1*fs1)*cos(11*pi*fd/3600) 
  710.             fy1 = fr*(1+f1*fs1)*sin(11*pi*fd/3600) 
  711.             aline fx1 fy1 
  712.             fx1 = fr*(1+f1*fs2)*cos(12*pi*fd/3600) 
  713.             fy1 = fr*(1+f1*fs2)*sin(12*pi*fd/3600) 
  714.             aline fx1 fy1 
  715.             fx1 = fr*(1+f1*fs3)*cos(13*pi*fd/3600) 
  716.             fy1 = fr*(1+f1*fs3)*sin(13*pi*fd/3600) 
  717.             aline fx1 fy1 
  718.             fx1 = fr*(1+f1*fs4)*cos(14*pi*fd/3600) 
  719.             fy1 = fr*(1+f1*fs4)*sin(14*pi*fd/3600) 
  720.             aline fx1 fy1 
  721.             fx1 = fr*(1+f1*fs5)*cos(15*pi*fd/3600) 
  722.             fy1 = fr*(1+f1*fs5)*sin(15*pi*fd/3600) 
  723.             aline fx1 fy1 
  724.             fx1 = fr*(1+f1*fs6)*cos(16*pi*fd/3600) 
  725.             fy1 = fr*(1+f1*fs6)*sin(16*pi*fd/3600) 
  726.             aline fx1 fy1 
  727.             fx1 = fr*(1+f1*fs7)*cos(17*pi*fd/3600) 
  728.             fy1 = fr*(1+f1*fs7)*sin(17*pi*fd/3600) 
  729.             aline fx1 fy1 
  730.             fx1 = fr*(1+f1*fs8)*cos(18*pi*fd/3600) 
  731.             fy1 = fr*(1+f1*fs8)*sin(18*pi*fd/3600) 
  732.             aline fx1 fy1 
  733.             fx1 = fr*(1+f1*fs9)*cos(19*pi*fd/3600) 
  734.             fy1 = fr*(1+f1*fs9)*sin(19*pi*fd/3600) 
  735.             aline fx1 fy1 
  736.             fx1 = fr*(1+f1)*cos(pi*fd/180) 
  737.             fy1 = fr*(1+f1)*sin(pi*fd/180) 
  738.             aline fx1 fy1 
  739.             end rotate 
  740.             f2 = f2-fd 
  741.             f1 = -f1 
  742.         end if 
  743.         begin rotate fa+fl 
  744.             fx1 = fr*(1-f1)*cos(pi*f2/180) 
  745.             fy1 = fr*(1-f1)*sin(pi*f2/180) 
  746.             amove fx1 -fy1 
  747.             fx1 = fr*(1-f1*fs9)*cos(9*pi*f2/1800) 
  748.             fy1 = fr*(1-f1*fs9)*sin(9*pi*f2/1800) 
  749.             aline fx1 -fy1 
  750.             fx1 = fr*(1-f1*fs8)*cos(8*pi*f2/1800) 
  751.             fy1 = fr*(1-f1*fs8)*sin(8*pi*f2/1800) 
  752.             aline fx1 -fy1 
  753.             fx1 = fr*(1-f1*fs7)*cos(7*pi*f2/1800) 
  754.             fy1 = fr*(1-f1*fs7)*sin(7*pi*f2/1800) 
  755.             aline fx1 -fy1 
  756.             fx1 = fr*(1-f1*fs6)*cos(6*pi*f2/1800) 
  757.             fy1 = fr*(1-f1*fs6)*sin(6*pi*f2/1800) 
  758.             aline fx1 -fy1 
  759.             fx1 = fr*(1-f1*fs5)*cos(5*pi*f2/1800) 
  760.             fy1 = fr*(1-f1*fs5)*sin(5*pi*f2/1800) 
  761.             aline fx1 -fy1 
  762.             fx1 = fr*(1-f1*fs4)*cos(4*pi*f2/1800) 
  763.             fy1 = fr*(1-f1*fs4)*sin(4*pi*f2/1800) 
  764.             aline fx1 -fy1 
  765.             fx1 = fr*(1-f1*fs3)*cos(3*pi*f2/1800) 
  766.             fy1 = fr*(1-f1*fs3)*sin(3*pi*f2/1800) 
  767.             aline fx1 -fy1 
  768.             fx1 = fr*(1-f1*fs2)*cos(2*pi*f2/1800) 
  769.             fy1 = fr*(1-f1*fs2)*sin(2*pi*f2/1800) 
  770.             aline fx1 -fy1 
  771.             fx1 = fr*(1-f1*fs1)*cos(1*pi*f2/1800) 
  772.             fy1 = fr*(1-f1*fs1)*sin(1*pi*f2/1800) 
  773.             aline fx1 -fy1 
  774.             aline fr 0 
  775.         end rotate 
  776.             else 
  777.         begin rotate fa 
  778.             amove fr 0 
  779.             fx1 = fr*(1-f1*fs1)*cos(1*pi*f2/1800) 
  780.             fy1 = fr*(1-f1*fs1)*sin(1*pi*f2/1800) 
  781.             aline fx1 -fy1 
  782.                 fx1 = fr*(1-f1*fs2)*cos(2*pi*f2/1800) 
  783.             fy1 = fr*(1-f1*fs2)*sin(2*pi*f2/1800) 
  784.             aline fx1 -fy1 
  785.                 fx1 = fr*(1-f1*fs3)*cos(3*pi*f2/1800) 
  786.             fy1 = fr*(1-f1*fs3)*sin(3*pi*f2/1800) 
  787.             aline fx1 -fy1 
  788.                 fx1 = fr*(1-f1*fs4)*cos(4*pi*f2/1800) 
  789.             fy1 = fr*(1-f1*fs4)*sin(4*pi*f2/1800) 
  790.             aline fx1 -fy1 
  791.                 fx1 = fr*(1-f1*fs5)*cos(5*pi*f2/1800) 
  792.             fy1 = fr*(1-f1*fs5)*sin(5*pi*f2/1800) 
  793.             aline fx1 -fy1 
  794.                 fx1 = fr*(1-f1*fs6)*cos(6*pi*f2/1800) 
  795.             fy1 = fr*(1-f1*fs6)*sin(6*pi*f2/1800) 
  796.             aline fx1 -fy1 
  797.                 fx1 = fr*(1-f1*fs7)*cos(7*pi*f2/1800) 
  798.             fy1 = fr*(1-f1*fs7)*sin(7*pi*f2/1800) 
  799.             aline fx1 -fy1 
  800.                 fx1 = fr*(1-f1*fs8)*cos(8*pi*f2/1800) 
  801.             fy1 = fr*(1-f1*fs8)*sin(8*pi*f2/1800) 
  802.             aline fx1 -fy1 
  803.                 fx1 = fr*(1-f1*fs9)*cos(9*pi*f2/1800) 
  804.             fy1 = fr*(1-f1*fs9)*sin(9*pi*f2/1800) 
  805.             aline fx1 -fy1 
  806.                 fx1 = fr*(1-f1)*cos(pi*f2/180) 
  807.             fy1 = fr*(1-f1)*sin(pi*f2/180) 
  808.             aline fx1 -fy1 
  809.         end rotate 
  810.         if fn>2.5 then 
  811.             for fi = f2 to fm step 2*fd 
  812.                 begin rotate fa-fi 
  813.                 amove fr*(1-f1) 0 
  814.                 fx1 = fr*(1-f1*fs9)*cos(1*pi*fd/3600) 
  815.                 fy1 = fr*(1-f1*fs9)*sin(1*pi*fd/3600) 
  816.                 aline fx1 -fy1 
  817.                 fx1 = fr*(1-f1*fs8)*cos(2*pi*fd/3600) 
  818.                 fy1 = fr*(1-f1*fs8)*sin(2*pi*fd/3600) 
  819.                 aline fx1 -fy1 
  820.                 fx1 = fr*(1-f1*fs7)*cos(3*pi*fd/3600) 
  821.                 fy1 = fr*(1-f1*fs7)*sin(3*pi*fd/3600) 
  822.                 aline fx1 -fy1 
  823.                 fx1 = fr*(1-f1*fs6)*cos(4*pi*fd/3600) 
  824.                 fy1 = fr*(1-f1*fs6)*sin(4*pi*fd/3600) 
  825.                 aline fx1 -fy1 
  826.                 fx1 = fr*(1-f1*fs5)*cos(5*pi*fd/3600) 
  827.                 fy1 = fr*(1-f1*fs5)*sin(5*pi*fd/3600) 
  828.                 aline fx1 -fy1 
  829.                 fx1 = fr*(1-f1*fs4)*cos(6*pi*fd/3600) 
  830.                 fy1 = fr*(1-f1*fs4)*sin(6*pi*fd/3600) 
  831.                 aline fx1 -fy1 
  832.                 fx1 = fr*(1-f1*fs3)*cos(7*pi*fd/3600) 
  833.                 fy1 = fr*(1-f1*fs3)*sin(7*pi*fd/3600) 
  834.                 aline fx1 -fy1 
  835.                 fx1 = fr*(1-f1*fs2)*cos(8*pi*fd/3600) 
  836.                 fy1 = fr*(1-f1*fs2)*sin(8*pi*fd/3600) 
  837.                 aline fx1 -fy1 
  838.                 fx1 = fr*(1-f1*fs1)*cos(9*pi*fd/3600) 
  839.                 fy1 = fr*(1-f1*fs1)*sin(9*pi*fd/3600) 
  840.                 aline fx1 -fy1 
  841.                 fx1 = fr*cos(pi*fd/360) 
  842.                 fy1 = fr*sin(pi*fd/360) 
  843.                 aline fx1 -fy1 
  844.                 fx1 = fr*(1+f1*fs1)*cos(11*pi*fd/3600) 
  845.                 fy1 = fr*(1+f1*fs1)*sin(11*pi*fd/3600) 
  846.                 aline fx1 -fy1 
  847.                 fx1 = fr*(1+f1*fs2)*cos(12*pi*fd/3600) 
  848.                 fy1 = fr*(1+f1*fs2)*sin(12*pi*fd/3600) 
  849.                 aline fx1 -fy1 
  850.                 fx1 = fr*(1+f1*fs3)*cos(13*pi*fd/3600) 
  851.                 fy1 = fr*(1+f1*fs3)*sin(13*pi*fd/3600) 
  852.                 aline fx1 -fy1 
  853.                 fx1 = fr*(1+f1*fs4)*cos(14*pi*fd/3600) 
  854.                 fy1 = fr*(1+f1*fs4)*sin(14*pi*fd/3600) 
  855.                 aline fx1 -fy1 
  856.                 fx1 = fr*(1+f1*fs5)*cos(15*pi*fd/3600) 
  857.                 fy1 = fr*(1+f1*fs5)*sin(15*pi*fd/3600) 
  858.                 aline fx1 -fy1 
  859.                 fx1 = fr*(1+f1*fs6)*cos(16*pi*fd/3600) 
  860.                 fy1 = fr*(1+f1*fs6)*sin(16*pi*fd/3600) 
  861.                 aline fx1 -fy1 
  862.                 fx1 = fr*(1+f1*fs7)*cos(17*pi*fd/3600) 
  863.                 fy1 = fr*(1+f1*fs7)*sin(17*pi*fd/3600) 
  864.                 aline fx1 -fy1 
  865.                 fx1 = fr*(1+f1*fs8)*cos(18*pi*fd/3600) 
  866.                 fy1 = fr*(1+f1*fs8)*sin(18*pi*fd/3600) 
  867.                 aline fx1 -fy1 
  868.                 fx1 = fr*(1+f1*fs9)*cos(19*pi*fd/3600) 
  869.                 fy1 = fr*(1+f1*fs9)*sin(19*pi*fd/3600) 
  870.                 aline fx1 -fy1 
  871.                 fx1 = fr*(1+f1)*cos(pi*fd/180) 
  872.                 fy1 = fr*(1+f1)*sin(pi*fd/180) 
  873.                 aline fx1 -fy1 
  874.                 fx1 = fr*(1+f1*fs9)*cos(21*pi*fd/3600) 
  875.                 fy1 = fr*(1+f1*fs9)*sin(21*pi*fd/3600) 
  876.                 aline fx1 -fy1 
  877.                 fx1 = fr*(1+f1*fs8)*cos(22*pi*fd/3600) 
  878.                 fy1 = fr*(1+f1*fs8)*sin(22*pi*fd/3600) 
  879.                 aline fx1 -fy1 
  880.                 fx1 = fr*(1+f1*fs7)*cos(23*pi*fd/3600) 
  881.                 fy1 = fr*(1+f1*fs7)*sin(23*pi*fd/3600) 
  882.                 aline fx1 -fy1 
  883.                 fx1 = fr*(1+f1*fs6)*cos(24*pi*fd/3600) 
  884.                 fy1 = fr*(1+f1*fs6)*sin(24*pi*fd/3600) 
  885.                 aline fx1 -fy1 
  886.                 fx1 = fr*(1+f1*fs5)*cos(25*pi*fd/3600) 
  887.                 fy1 = fr*(1+f1*fs5)*sin(25*pi*fd/3600) 
  888.                 aline fx1 -fy1 
  889.                 fx1 = fr*(1+f1*fs4)*cos(26*pi*fd/3600) 
  890.                 fy1 = fr*(1+f1*fs4)*sin(26*pi*fd/3600) 
  891.                 aline fx1 -fy1 
  892.                 fx1 = fr*(1+f1*fs3)*cos(27*pi*fd/3600) 
  893.                 fy1 = fr*(1+f1*fs3)*sin(27*pi*fd/3600) 
  894.                 aline fx1 -fy1 
  895.                 fx1 = fr*(1+f1*fs2)*cos(28*pi*fd/3600) 
  896.                 fy1 = fr*(1+f1*fs2)*sin(28*pi*fd/3600) 
  897.                 aline fx1 -fy1 
  898.                 fx1 = fr*(1+f1*fs1)*cos(29*pi*fd/3600) 
  899.                 fy1 = fr*(1+f1*fs1)*sin(29*pi*fd/3600) 
  900.                 aline fx1 -fy1 
  901.                 fx1 = fr*cos(pi*fd/120) 
  902.                 fy1 = fr*sin(pi*fd/120) 
  903.                 aline fx1 -fy1 
  904.                 fx1 = fr*(1-f1*fs1)*cos(31*pi*fd/3600) 
  905.                 fy1 = fr*(1-f1*fs1)*sin(31*pi*fd/3600) 
  906.                 aline fx1 -fy1 
  907.                 fx1 = fr*(1-f1*fs2)*cos(32*pi*fd/3600) 
  908.                 fy1 = fr*(1-f1*fs2)*sin(32*pi*fd/3600) 
  909.                 aline fx1 -fy1 
  910.                 fx1 = fr*(1-f1*fs3)*cos(33*pi*fd/3600) 
  911.                 fy1 = fr*(1-f1*fs3)*sin(33*pi*fd/3600) 
  912.                 aline fx1 -fy1 
  913.                 fx1 = fr*(1-f1*fs4)*cos(34*pi*fd/3600) 
  914.                 fy1 = fr*(1-f1*fs4)*sin(34*pi*fd/3600) 
  915.                 aline fx1 -fy1 
  916.                 fx1 = fr*(1-f1*fs5)*cos(35*pi*fd/3600) 
  917.                 fy1 = fr*(1-f1*fs5)*sin(35*pi*fd/3600) 
  918.                 aline fx1 -fy1 
  919.                 fx1 = fr*(1-f1*fs6)*cos(36*pi*fd/3600) 
  920.                 fy1 = fr*(1-f1*fs6)*sin(36*pi*fd/3600) 
  921.                 aline fx1 -fy1 
  922.                 fx1 = fr*(1-f1*fs7)*cos(37*pi*fd/3600) 
  923.                 fy1 = fr*(1-f1*fs7)*sin(37*pi*fd/3600) 
  924.                 aline fx1 -fy1 
  925.                 fx1 = fr*(1-f1*fs8)*cos(38*pi*fd/3600) 
  926.                 fy1 = fr*(1-f1*fs8)*sin(38*pi*fd/3600) 
  927.                 aline fx1 -fy1 
  928.                 fx1 = fr*(1-f1*fs9)*cos(39*pi*fd/3600) 
  929.                 fy1 = fr*(1-f1*fs9)*sin(39*pi*fd/3600) 
  930.                 aline fx1 -fy1 
  931.                 fx1 = fr*(1-f1)*cos(pi*fd/90) 
  932.                 fy1 = fr*(1-f1)*sin(pi*fd/90) 
  933.                 aline fx1 -fy1 
  934.             end rotate 
  935.             next fi 
  936.         end if 
  937.         f2 = fl-fm-fd 
  938.         if abs(fn-2*fix(0.5*fn))<0.1 then 
  939.             begin rotate fa-fm-fd 
  940.             amove fr*(1-f1) 0 
  941.             fx1 = fr*(1-f1*fs9)*cos(1*pi*fd/3600) 
  942.             fy1 = fr*(1-f1*fs9)*sin(1*pi*fd/3600) 
  943.             aline fx1 -fy1 
  944.             fx1 = fr*(1-f1*fs8)*cos(2*pi*fd/3600) 
  945.             fy1 = fr*(1-f1*fs8)*sin(2*pi*fd/3600) 
  946.             aline fx1 -fy1 
  947.             fx1 = fr*(1-f1*fs7)*cos(3*pi*fd/3600) 
  948.             fy1 = fr*(1-f1*fs7)*sin(3*pi*fd/3600) 
  949.             aline fx1 -fy1 
  950.             fx1 = fr*(1-f1*fs6)*cos(4*pi*fd/3600) 
  951.             fy1 = fr*(1-f1*fs6)*sin(4*pi*fd/3600) 
  952.             aline fx1 -fy1 
  953.             fx1 = fr*(1-f1*fs5)*cos(5*pi*fd/3600) 
  954.             fy1 = fr*(1-f1*fs5)*sin(5*pi*fd/3600) 
  955.             aline fx1 -fy1 
  956.             fx1 = fr*(1-f1*fs4)*cos(6*pi*fd/3600) 
  957.             fy1 = fr*(1-f1*fs4)*sin(6*pi*fd/3600) 
  958.             aline fx1 -fy1 
  959.             fx1 = fr*(1-f1*fs3)*cos(7*pi*fd/3600) 
  960.             fy1 = fr*(1-f1*fs3)*sin(7*pi*fd/3600) 
  961.             aline fx1 -fy1 
  962.             fx1 = fr*(1-f1*fs2)*cos(8*pi*fd/3600) 
  963.             fy1 = fr*(1-f1*fs2)*sin(8*pi*fd/3600) 
  964.             aline fx1 -fy1 
  965.             fx1 = fr*(1-f1*fs1)*cos(9*pi*fd/3600) 
  966.             fy1 = fr*(1-f1*fs1)*sin(9*pi*fd/3600) 
  967.             aline fx1 -fy1 
  968.             fx1 = fr*cos(pi*fd/360) 
  969.             fy1 = fr*sin(pi*fd/360) 
  970.             aline fx1 -fy1 
  971.             fx1 = fr*(1+f1*fs1)*cos(11*pi*fd/3600) 
  972.             fy1 = fr*(1+f1*fs1)*sin(11*pi*fd/3600) 
  973.             aline fx1 -fy1 
  974.             fx1 = fr*(1+f1*fs2)*cos(12*pi*fd/3600) 
  975.             fy1 = fr*(1+f1*fs2)*sin(12*pi*fd/3600) 
  976.             aline fx1 -fy1 
  977.             fx1 = fr*(1+f1*fs3)*cos(13*pi*fd/3600) 
  978.             fy1 = fr*(1+f1*fs3)*sin(13*pi*fd/3600) 
  979.             aline fx1 -fy1 
  980.             fx1 = fr*(1+f1*fs4)*cos(14*pi*fd/3600) 
  981.             fy1 = fr*(1+f1*fs4)*sin(14*pi*fd/3600) 
  982.             aline fx1 -fy1 
  983.             fx1 = fr*(1+f1*fs5)*cos(15*pi*fd/3600) 
  984.             fy1 = fr*(1+f1*fs5)*sin(15*pi*fd/3600) 
  985.             aline fx1 -fy1 
  986.             fx1 = fr*(1+f1*fs6)*cos(16*pi*fd/3600) 
  987.             fy1 = fr*(1+f1*fs6)*sin(16*pi*fd/3600) 
  988.             aline fx1 -fy1 
  989.             fx1 = fr*(1+f1*fs7)*cos(17*pi*fd/3600) 
  990.             fy1 = fr*(1+f1*fs7)*sin(17*pi*fd/3600) 
  991.             aline fx1 -fy1 
  992.             fx1 = fr*(1+f1*fs8)*cos(18*pi*fd/3600) 
  993.             fy1 = fr*(1+f1*fs8)*sin(18*pi*fd/3600) 
  994.             aline fx1 -fy1 
  995.             fx1 = fr*(1+f1*fs9)*cos(19*pi*fd/3600) 
  996.             fy1 = fr*(1+f1*fs9)*sin(19*pi*fd/3600) 
  997.             aline fx1 -fy1 
  998.             fx1 = fr*(1+f1)*cos(pi*fd/180) 
  999.             fy1 = fr*(1+f1)*sin(pi*fd/180) 
  1000.             aline fx1 -fy1 
  1001.             end rotate 
  1002.             f2 = f2-fd 
  1003.             f1 = -f1 
  1004.         end if 
  1005.         begin rotate fa-fl 
  1006.             fx1 = fr*(1-f1)*cos(pi*f2/180) 
  1007.             fy1 = fr*(1-f1)*sin(pi*f2/180) 
  1008.             amove fx1 fy1 
  1009.             fx1 = fr*(1-f1*fs9)*cos(9*pi*f2/1800) 
  1010.             fy1 = fr*(1-f1*fs9)*sin(9*pi*f2/1800) 
  1011.             aline fx1 fy1 
  1012.             fx1 = fr*(1-f1*fs8)*cos(8*pi*f2/1800) 
  1013.             fy1 = fr*(1-f1*fs8)*sin(8*pi*f2/1800) 
  1014.             aline fx1 fy1 
  1015.             fx1 = fr*(1-f1*fs7)*cos(7*pi*f2/1800) 
  1016.             fy1 = fr*(1-f1*fs7)*sin(7*pi*f2/1800) 
  1017.             aline fx1 fy1 
  1018.             fx1 = fr*(1-f1*fs6)*cos(6*pi*f2/1800) 
  1019.             fy1 = fr*(1-f1*fs6)*sin(6*pi*f2/1800) 
  1020.             aline fx1 fy1 
  1021.             fx1 = fr*(1-f1*fs5)*cos(5*pi*f2/1800) 
  1022.             fy1 = fr*(1-f1*fs5)*sin(5*pi*f2/1800) 
  1023.             aline fx1 fy1 
  1024.             fx1 = fr*(1-f1*fs4)*cos(4*pi*f2/1800) 
  1025.             fy1 = fr*(1-f1*fs4)*sin(4*pi*f2/1800) 
  1026.             aline fx1 fy1 
  1027.             fx1 = fr*(1-f1*fs3)*cos(3*pi*f2/1800) 
  1028.             fy1 = fr*(1-f1*fs3)*sin(3*pi*f2/1800) 
  1029.             aline fx1 fy1 
  1030.             fx1 = fr*(1-f1*fs2)*cos(2*pi*f2/1800) 
  1031.             fy1 = fr*(1-f1*fs2)*sin(2*pi*f2/1800) 
  1032.             aline fx1 fy1 
  1033.             fx1 = fr*(1-f1*fs1)*cos(1*pi*f2/1800) 
  1034.             fy1 = fr*(1-f1*fs1)*sin(1*pi*f2/1800) 
  1035.             aline fx1 fy1 
  1036.             aline fr 0 
  1037.         end rotate 
  1038.             end if 
  1039.         end origin 
  1040.     amove ffx ffy 
  1041.     end if 
  1042. end sub 
  1043.  
  1044. GluonW = 1  ! Winding coefficient, should be >1/pi 
  1045.  
  1046. ! Gluon line from the current point to x,y 
  1047. sub Gluon x y 
  1048.     @feyn x y 
  1049.     if fl>0 then 
  1050.         fn = fix(0.5*fl/PhotonL-GluonW) 
  1051.         if fn<1 then 
  1052.             fn = 1 
  1053.         end if 
  1054.         fd = fl/(2*fn+1+2*GluonW) 
  1055.         fm = fd*(2*fn-1) 
  1056.         begin origin 
  1057.             begin rotate fa 
  1058.                 for fi = 0 to fm step 2*fd 
  1059.                     begin translate fi 0 
  1060.             amove 0 0 
  1061.             aline (0.05+GluonW*(1-fs9))*fd  fs1*PhotonA 
  1062.             aline (0.1+GluonW*(1-fs8))*fd   fs2*PhotonA 
  1063.             aline (0.15+GluonW*(1-fs7))*fd  fs3*PhotonA 
  1064.             aline (0.2+GluonW*(1-fs6))*fd   fs4*PhotonA 
  1065.             aline (0.25+GluonW*(1-fs5))*fd  fs5*PhotonA 
  1066.             aline (0.3+GluonW*(1-fs4))*fd   fs6*PhotonA 
  1067.             aline (0.35+GluonW*(1-fs3))*fd  fs7*PhotonA 
  1068.             aline (0.4+GluonW*(1-fs2))*fd   fs8*PhotonA 
  1069.             aline (0.45+GluonW*(1-fs1))*fd  fs9*PhotonA 
  1070.             aline (0.5+GluonW)*fd           PhotonA 
  1071.             aline (0.55+GluonW*(1+fs1))*fd  fs9*PhotonA 
  1072.             aline (0.6+GluonW*(1+fs2))*fd   fs8*PhotonA 
  1073.             aline (0.65+GluonW*(1+fs3))*fd  fs7*PhotonA 
  1074.             aline (0.7+GluonW*(1+fs4))*fd   fs6*PhotonA 
  1075.             aline (0.75+GluonW*(1+fs5))*fd  fs5*PhotonA 
  1076.             aline (0.8+GluonW*(1+fs6))*fd   fs4*PhotonA 
  1077.             aline (0.85+GluonW*(1+fs7))*fd  fs3*PhotonA 
  1078.             aline (0.9+GluonW*(1+fs8))*fd   fs2*PhotonA 
  1079.             aline (0.95+GluonW*(1+fs9))*fd  fs1*PhotonA 
  1080.             aline (1+2*GluonW)*fd           0 
  1081.             aline (1.05+GluonW*(1+fs9))*fd -fs1*PhotonA 
  1082.             aline (1.1+GluonW*(1+fs8))*fd  -fs2*PhotonA 
  1083.             aline (1.15+GluonW*(1+fs7))*fd -fs3*PhotonA 
  1084.             aline (1.2+GluonW*(1+fs6))*fd  -fs4*PhotonA 
  1085.             aline (1.25+GluonW*(1+fs5))*fd -fs5*PhotonA 
  1086.             aline (1.3+GluonW*(1+fs4))*fd  -fs6*PhotonA 
  1087.             aline (1.35+GluonW*(1+fs3))*fd -fs7*PhotonA 
  1088.             aline (1.4+GluonW*(1+fs2))*fd  -fs8*PhotonA 
  1089.             aline (1.45+GluonW*(1+fs1))*fd -fs9*PhotonA 
  1090.             aline (1.5+GluonW)*fd          -PhotonA 
  1091.             aline (1.55+GluonW*(1-fs1))*fd -fs9*PhotonA 
  1092.             aline (1.6+GluonW*(1-fs2))*fd  -fs8*PhotonA 
  1093.             aline (1.65+GluonW*(1-fs3))*fd -fs7*PhotonA 
  1094.             aline (1.7+GluonW*(1-fs4))*fd  -fs6*PhotonA 
  1095.             aline (1.75+GluonW*(1-fs5))*fd -fs5*PhotonA 
  1096.             aline (1.8+GluonW*(1-fs6))*fd  -fs4*PhotonA 
  1097.             aline (1.85+GluonW*(1-fs7))*fd -fs3*PhotonA 
  1098.             aline (1.9+GluonW*(1-fs8))*fd  -fs2*PhotonA 
  1099.             aline (1.95+GluonW*(1-fs9))*fd -fs1*PhotonA 
  1100.             aline 2*fd                      0 
  1101.                     end translate 
  1102.                 next fi 
  1103.         begin translate 2*fn*fd 0 
  1104.             amove 0 0 
  1105.             aline (0.05+GluonW*(1-fs9))*fd fs1*PhotonA 
  1106.             aline (0.1+GluonW*(1-fs8))*fd  fs2*PhotonA 
  1107.             aline (0.15+GluonW*(1-fs7))*fd fs3*PhotonA 
  1108.             aline (0.2+GluonW*(1-fs6))*fd  fs4*PhotonA 
  1109.             aline (0.25+GluonW*(1-fs5))*fd fs5*PhotonA 
  1110.             aline (0.3+GluonW*(1-fs4))*fd  fs6*PhotonA 
  1111.             aline (0.35+GluonW*(1-fs3))*fd fs7*PhotonA 
  1112.             aline (0.4+GluonW*(1-fs2))*fd  fs8*PhotonA 
  1113.             aline (0.45+GluonW*(1-fs1))*fd fs9*PhotonA 
  1114.             aline (0.5+GluonW)*fd          PhotonA 
  1115.             aline (0.55+GluonW*(1+fs1))*fd fs9*PhotonA 
  1116.             aline (0.6+GluonW*(1+fs2))*fd  fs8*PhotonA 
  1117.             aline (0.65+GluonW*(1+fs3))*fd fs7*PhotonA 
  1118.             aline (0.7+GluonW*(1+fs4))*fd  fs6*PhotonA 
  1119.             aline (0.75+GluonW*(1+fs5))*fd fs5*PhotonA 
  1120.             aline (0.8+GluonW*(1+fs6))*fd  fs4*PhotonA 
  1121.             aline (0.85+GluonW*(1+fs7))*fd fs3*PhotonA 
  1122.             aline (0.9+GluonW*(1+fs8))*fd  fs2*PhotonA 
  1123.             aline (0.95+GluonW*(1+fs9))*fd fs1*PhotonA 
  1124.             aline (1+2*GluonW)*fd          0 
  1125.         end translate 
  1126.             end rotate 
  1127.         end origin 
  1128.     amove ffx ffy 
  1129.     end if 
  1130. end sub 
  1131.  
  1132. ! Gluon line from the current point via xx,yy to x,y 
  1133. sub Gluon2 xx yy x y 
  1134.     @feyn2 xx yy x y 
  1135.     if fr<PhotonA+fe then 
  1136.         @Gluon x y 
  1137.     else 
  1138.         fn = fix(fr*fl*pi/360/PhotonL-GluonW) 
  1139.         if fn<1 then 
  1140.             fn = 1 
  1141.         end if 
  1142.         fd = fl/(2*fn+1+2*GluonW) 
  1143.     fm = fd*(2*fn-1) 
  1144.         amove fx fy 
  1145.         begin origin 
  1146.             if fc>0 then 
  1147.         for fi = 0 to fm step 2*fd 
  1148.             begin rotate fa+fi 
  1149.             amove fr 0 
  1150.             f1 = fr-PhotonA*fs1 
  1151.             f2 = (0.05+GluonW*(1-fs9))*fd/180*pi 
  1152.             aline f1*cos(f2) f1*sin(f2) 
  1153.             f1 = fr-PhotonA*fs2 
  1154.             f2 = (0.1+GluonW*(1-fs8))*fd/180*pi 
  1155.             aline f1*cos(f2) f1*sin(f2) 
  1156.             f1 = fr-PhotonA*fs3 
  1157.             f2 = (0.15+GluonW*(1-fs7))*fd/180*pi 
  1158.             aline f1*cos(f2) f1*sin(f2) 
  1159.             f1 = fr-PhotonA*fs4 
  1160.             f2 = (0.2+GluonW*(1-fs6))*fd/180*pi 
  1161.             aline f1*cos(f2) f1*sin(f2) 
  1162.             f1 = fr-PhotonA*fs5 
  1163.             f2 = (0.25+GluonW*(1-fs5))*fd/180*pi 
  1164.             aline f1*cos(f2) f1*sin(f2) 
  1165.             f1 = fr-PhotonA*fs6 
  1166.             f2 = (0.3+GluonW*(1-fs4))*fd/180*pi 
  1167.             aline f1*cos(f2) f1*sin(f2) 
  1168.             f1 = fr-PhotonA*fs7 
  1169.             f2 = (0.35+GluonW*(1-fs3))*fd/180*pi 
  1170.             aline f1*cos(f2) f1*sin(f2) 
  1171.             f1 = fr-PhotonA*fs8 
  1172.             f2 = (0.4+GluonW*(1-fs2))*fd/180*pi 
  1173.             aline f1*cos(f2) f1*sin(f2) 
  1174.             f1 = fr-PhotonA*fs9 
  1175.             f2 = (0.45+GluonW*(1-fs1))*fd/180*pi 
  1176.             aline f1*cos(f2) f1*sin(f2) 
  1177.             f1 = fr-PhotonA 
  1178.             f2 = (0.5+GluonW)*fd/180*pi 
  1179.             aline f1*cos(f2) f1*sin(f2) 
  1180.             f1 = fr-PhotonA*fs9 
  1181.             f2 = (0.55+GluonW*(1+fs1))*fd/180*pi 
  1182.             aline f1*cos(f2) f1*sin(f2) 
  1183.             f1 = fr-PhotonA*fs8 
  1184.             f2 = (0.6+GluonW*(1+fs2))*fd/180*pi 
  1185.             aline f1*cos(f2) f1*sin(f2) 
  1186.             f1 = fr-PhotonA*fs7 
  1187.             f2 = (0.65+GluonW*(1+fs3))*fd/180*pi 
  1188.             aline f1*cos(f2) f1*sin(f2) 
  1189.             f1 = fr-PhotonA*fs6 
  1190.             f2 = (0.7+GluonW*(1+fs4))*fd/180*pi 
  1191.             aline f1*cos(f2) f1*sin(f2) 
  1192.             f1 = fr-PhotonA*fs5 
  1193.             f2 = (0.75+GluonW*(1+fs5))*fd/180*pi 
  1194.             aline f1*cos(f2) f1*sin(f2) 
  1195.             f1 = fr-PhotonA*fs4 
  1196.             f2 = (0.8+GluonW*(1+fs6))*fd/180*pi 
  1197.             aline f1*cos(f2) f1*sin(f2) 
  1198.             f1 = fr-PhotonA*fs3 
  1199.             f2 = (0.85+GluonW*(1+fs7))*fd/180*pi 
  1200.             aline f1*cos(f2) f1*sin(f2) 
  1201.             f1 = fr-PhotonA*fs2 
  1202.             f2 = (0.9+GluonW*(1+fs8))*fd/180*pi 
  1203.             aline f1*cos(f2) f1*sin(f2) 
  1204.             f1 = fr-PhotonA*fs1 
  1205.             f2 = (0.95+GluonW*(1+fs9))*fd/180*pi 
  1206.             aline f1*cos(f2) f1*sin(f2) 
  1207.             f1 = fr 
  1208.             f2 = (1+2*GluonW)*fd/180*pi 
  1209.             aline f1*cos(f2) f1*sin(f2) 
  1210.             f1 = fr+PhotonA*fs1 
  1211.             f2 = (1.05+GluonW*(1+fs9))*fd/180*pi 
  1212.             aline f1*cos(f2) f1*sin(f2) 
  1213.             f1 = fr+PhotonA*fs2 
  1214.             f2 = (1.1+GluonW*(1+fs8))*fd/180*pi 
  1215.             aline f1*cos(f2) f1*sin(f2) 
  1216.             f1 = fr+PhotonA*fs3 
  1217.             f2 = (1.15+GluonW*(1+fs7))*fd/180*pi 
  1218.             aline f1*cos(f2) f1*sin(f2) 
  1219.             f1 = fr+PhotonA*fs4 
  1220.             f2 = (1.2+GluonW*(1+fs6))*fd/180*pi 
  1221.             aline f1*cos(f2) f1*sin(f2) 
  1222.             f1 = fr+PhotonA*fs5 
  1223.             f2 = (1.25+GluonW*(1+fs5))*fd/180*pi 
  1224.             aline f1*cos(f2) f1*sin(f2) 
  1225.             f1 = fr+PhotonA*fs6 
  1226.             f2 = (1.3+GluonW*(1+fs4))*fd/180*pi 
  1227.             aline f1*cos(f2) f1*sin(f2) 
  1228.             f1 = fr+PhotonA*fs7 
  1229.             f2 = (1.35+GluonW*(1+fs3))*fd/180*pi 
  1230.             aline f1*cos(f2) f1*sin(f2) 
  1231.             f1 = fr+PhotonA*fs8 
  1232.             f2 = (1.4+GluonW*(1+fs2))*fd/180*pi 
  1233.             aline f1*cos(f2) f1*sin(f2) 
  1234.             f1 = fr+PhotonA*fs9 
  1235.             f2 = (1.45+GluonW*(1+fs1))*fd/180*pi 
  1236.             aline f1*cos(f2) f1*sin(f2) 
  1237.             f1 = fr+PhotonA 
  1238.             f2 = (1.5+GluonW)*fd/180*pi 
  1239.             aline f1*cos(f2) f1*sin(f2) 
  1240.             f1 = fr+PhotonA*fs9 
  1241.             f2 = (1.55+GluonW*(1-fs1))*fd/180*pi 
  1242.             aline f1*cos(f2) f1*sin(f2) 
  1243.             f1 = fr+PhotonA*fs8 
  1244.             f2 = (1.6+GluonW*(1-fs2))*fd/180*pi 
  1245.             aline f1*cos(f2) f1*sin(f2) 
  1246.             f1 = fr+PhotonA*fs7 
  1247.             f2 = (1.65+GluonW*(1-fs3))*fd/180*pi 
  1248.             aline f1*cos(f2) f1*sin(f2) 
  1249.             f1 = fr+PhotonA*fs6 
  1250.             f2 = (1.7+GluonW*(1-fs4))*fd/180*pi 
  1251.             aline f1*cos(f2) f1*sin(f2) 
  1252.             f1 = fr+PhotonA*fs5 
  1253.             f2 = (1.75+GluonW*(1-fs5))*fd/180*pi 
  1254.             aline f1*cos(f2) f1*sin(f2) 
  1255.             f1 = fr+PhotonA*fs4 
  1256.             f2 = (1.8+GluonW*(1-fs6))*fd/180*pi 
  1257.             aline f1*cos(f2) f1*sin(f2) 
  1258.             f1 = fr+PhotonA*fs3 
  1259.             f2 = (1.85+GluonW*(1-fs7))*fd/180*pi 
  1260.             aline f1*cos(f2) f1*sin(f2) 
  1261.             f1 = fr+PhotonA*fs2 
  1262.             f2 = (1.9+GluonW*(1-fs8))*fd/180*pi 
  1263.             aline f1*cos(f2) f1*sin(f2) 
  1264.             f1 = fr+PhotonA*fs1 
  1265.             f2 = (1.95+GluonW*(1-fs9))*fd/180*pi 
  1266.             aline f1*cos(f2) f1*sin(f2) 
  1267.             f1 = fr 
  1268.             f2 = 2*fd/180*pi 
  1269.             aline f1*cos(f2) f1*sin(f2) 
  1270.             end rotate 
  1271.         next fi 
  1272.         begin rotate fa+2*fn*fd 
  1273.             amove fr 0 
  1274.             f1 = fr-PhotonA*fs1 
  1275.             f2 = (0.05+GluonW*(1-fs9))*fd/180*pi 
  1276.             aline f1*cos(f2) f1*sin(f2) 
  1277.             f1 = fr-PhotonA*fs2 
  1278.             f2 = (0.1+GluonW*(1-fs8))*fd/180*pi 
  1279.             aline f1*cos(f2) f1*sin(f2) 
  1280.             f1 = fr-PhotonA*fs3 
  1281.             f2 = (0.15+GluonW*(1-fs7))*fd/180*pi 
  1282.             aline f1*cos(f2) f1*sin(f2) 
  1283.             f1 = fr-PhotonA*fs4 
  1284.             f2 = (0.2+GluonW*(1-fs6))*fd/180*pi 
  1285.             aline f1*cos(f2) f1*sin(f2) 
  1286.             f1 = fr-PhotonA*fs5 
  1287.             f2 = (0.25+GluonW*(1-fs5))*fd/180*pi 
  1288.             aline f1*cos(f2) f1*sin(f2) 
  1289.             f1 = fr-PhotonA*fs6 
  1290.             f2 = (0.3+GluonW*(1-fs4))*fd/180*pi 
  1291.             aline f1*cos(f2) f1*sin(f2) 
  1292.             f1 = fr-PhotonA*fs7 
  1293.             f2 = (0.35+GluonW*(1-fs3))*fd/180*pi 
  1294.             aline f1*cos(f2) f1*sin(f2) 
  1295.             f1 = fr-PhotonA*fs8 
  1296.             f2 = (0.4+GluonW*(1-fs2))*fd/180*pi 
  1297.             aline f1*cos(f2) f1*sin(f2) 
  1298.             f1 = fr-PhotonA*fs9 
  1299.             f2 = (0.45+GluonW*(1-fs1))*fd/180*pi 
  1300.             aline f1*cos(f2) f1*sin(f2) 
  1301.             f1 = fr-PhotonA 
  1302.             f2 = (0.5+GluonW)*fd/180*pi 
  1303.             aline f1*cos(f2) f1*sin(f2) 
  1304.             f1 = fr-PhotonA*fs9 
  1305.             f2 = (0.55+GluonW*(1+fs1))*fd/180*pi 
  1306.             aline f1*cos(f2) f1*sin(f2) 
  1307.             f1 = fr-PhotonA*fs8 
  1308.             f2 = (0.6+GluonW*(1+fs2))*fd/180*pi 
  1309.             aline f1*cos(f2) f1*sin(f2) 
  1310.             f1 = fr-PhotonA*fs7 
  1311.             f2 = (0.65+GluonW*(1+fs3))*fd/180*pi 
  1312.             aline f1*cos(f2) f1*sin(f2) 
  1313.             f1 = fr-PhotonA*fs6 
  1314.             f2 = (0.7+GluonW*(1+fs4))*fd/180*pi 
  1315.             aline f1*cos(f2) f1*sin(f2) 
  1316.             f1 = fr-PhotonA*fs5 
  1317.             f2 = (0.75+GluonW*(1+fs5))*fd/180*pi 
  1318.             aline f1*cos(f2) f1*sin(f2) 
  1319.             f1 = fr-PhotonA*fs4 
  1320.             f2 = (0.8+GluonW*(1+fs6))*fd/180*pi 
  1321.             aline f1*cos(f2) f1*sin(f2) 
  1322.             f1 = fr-PhotonA*fs3 
  1323.             f2 = (0.85+GluonW*(1+fs7))*fd/180*pi 
  1324.             aline f1*cos(f2) f1*sin(f2) 
  1325.             f1 = fr-PhotonA*fs2 
  1326.             f2 = (0.9+GluonW*(1+fs8))*fd/180*pi 
  1327.             aline f1*cos(f2) f1*sin(f2) 
  1328.             f1 = fr-PhotonA*fs1 
  1329.             f2 = (0.95+GluonW*(1+fs9))*fd/180*pi 
  1330.             aline f1*cos(f2) f1*sin(f2) 
  1331.             f1 = fr 
  1332.             f2 = (1+2*GluonW)*fd/180*pi 
  1333.             aline f1*cos(f2) f1*sin(f2) 
  1334.         end rotate 
  1335.             else 
  1336.         for fi = 0 to fm step 2*fd 
  1337.             begin rotate fa-fi 
  1338.             amove fr 0 
  1339.             f1 = fr+PhotonA*fs1 
  1340.             f2 = (0.05+GluonW*(1-fs9))*fd/180*pi 
  1341.             aline f1*cos(f2) -f1*sin(f2) 
  1342.             f1 = fr+PhotonA*fs2 
  1343.             f2 = (0.1+GluonW*(1-fs8))*fd/180*pi 
  1344.             aline f1*cos(f2) -f1*sin(f2) 
  1345.             f1 = fr+PhotonA*fs3 
  1346.             f2 = (0.15+GluonW*(1-fs7))*fd/180*pi 
  1347.             aline f1*cos(f2) -f1*sin(f2) 
  1348.             f1 = fr+PhotonA*fs4 
  1349.             f2 = (0.2+GluonW*(1-fs6))*fd/180*pi 
  1350.             aline f1*cos(f2) -f1*sin(f2) 
  1351.             f1 = fr+PhotonA*fs5 
  1352.             f2 = (0.25+GluonW*(1-fs5))*fd/180*pi 
  1353.             aline f1*cos(f2) -f1*sin(f2) 
  1354.             f1 = fr+PhotonA*fs6 
  1355.             f2 = (0.3+GluonW*(1-fs4))*fd/180*pi 
  1356.             aline f1*cos(f2) -f1*sin(f2) 
  1357.             f1 = fr+PhotonA*fs7 
  1358.             f2 = (0.35+GluonW*(1-fs3))*fd/180*pi 
  1359.             aline f1*cos(f2) -f1*sin(f2) 
  1360.             f1 = fr+PhotonA*fs8 
  1361.             f2 = (0.4+GluonW*(1-fs2))*fd/180*pi 
  1362.             aline f1*cos(f2) -f1*sin(f2) 
  1363.             f1 = fr+PhotonA*fs9 
  1364.             f2 = (0.45+GluonW*(1-fs1))*fd/180*pi 
  1365.             aline f1*cos(f2) -f1*sin(f2) 
  1366.             f1 = fr+PhotonA 
  1367.             f2 = (0.5+GluonW)*fd/180*pi 
  1368.             aline f1*cos(f2) -f1*sin(f2) 
  1369.             f1 = fr+PhotonA*fs9 
  1370.             f2 = (0.55+GluonW*(1+fs1))*fd/180*pi 
  1371.             aline f1*cos(f2) -f1*sin(f2) 
  1372.             f1 = fr+PhotonA*fs8 
  1373.             f2 = (0.6+GluonW*(1+fs2))*fd/180*pi 
  1374.             aline f1*cos(f2) -f1*sin(f2) 
  1375.             f1 = fr+PhotonA*fs7 
  1376.             f2 = (0.65+GluonW*(1+fs3))*fd/180*pi 
  1377.             aline f1*cos(f2) -f1*sin(f2) 
  1378.             f1 = fr+PhotonA*fs6 
  1379.             f2 = (0.7+GluonW*(1+fs4))*fd/180*pi 
  1380.             aline f1*cos(f2) -f1*sin(f2) 
  1381.             f1 = fr+PhotonA*fs5 
  1382.             f2 = (0.75+GluonW*(1+fs5))*fd/180*pi 
  1383.             aline f1*cos(f2) -f1*sin(f2) 
  1384.             f1 = fr+PhotonA*fs4 
  1385.             f2 = (0.8+GluonW*(1+fs6))*fd/180*pi 
  1386.             aline f1*cos(f2) -f1*sin(f2) 
  1387.             f1 = fr+PhotonA*fs3 
  1388.             f2 = (0.85+GluonW*(1+fs7))*fd/180*pi 
  1389.             aline f1*cos(f2) -f1*sin(f2) 
  1390.             f1 = fr+PhotonA*fs2 
  1391.             f2 = (0.9+GluonW*(1+fs8))*fd/180*pi 
  1392.             aline f1*cos(f2) -f1*sin(f2) 
  1393.             f1 = fr+PhotonA*fs1 
  1394.             f2 = (0.95+GluonW*(1+fs9))*fd/180*pi 
  1395.             aline f1*cos(f2) -f1*sin(f2) 
  1396.             f1 = fr 
  1397.             f2 = (1+2*GluonW)*fd/180*pi 
  1398.             aline f1*cos(f2) -f1*sin(f2) 
  1399.             f1 = fr-PhotonA*fs1 
  1400.             f2 = (1.05+GluonW*(1+fs9))*fd/180*pi 
  1401.             aline f1*cos(f2) -f1*sin(f2) 
  1402.             f1 = fr-PhotonA*fs2 
  1403.             f2 = (1.1+GluonW*(1+fs8))*fd/180*pi 
  1404.             aline f1*cos(f2) -f1*sin(f2) 
  1405.             f1 = fr-PhotonA*fs3 
  1406.             f2 = (1.15+GluonW*(1+fs7))*fd/180*pi 
  1407.             aline f1*cos(f2) -f1*sin(f2) 
  1408.             f1 = fr-PhotonA*fs4 
  1409.             f2 = (1.2+GluonW*(1+fs6))*fd/180*pi 
  1410.             aline f1*cos(f2) -f1*sin(f2) 
  1411.             f1 = fr-PhotonA*fs5 
  1412.             f2 = (1.25+GluonW*(1+fs5))*fd/180*pi 
  1413.             aline f1*cos(f2) -f1*sin(f2) 
  1414.             f1 = fr-PhotonA*fs6 
  1415.             f2 = (1.3+GluonW*(1+fs4))*fd/180*pi 
  1416.             aline f1*cos(f2) -f1*sin(f2) 
  1417.             f1 = fr-PhotonA*fs7 
  1418.             f2 = (1.35+GluonW*(1+fs3))*fd/180*pi 
  1419.             aline f1*cos(f2) -f1*sin(f2) 
  1420.             f1 = fr-PhotonA*fs8 
  1421.             f2 = (1.4+GluonW*(1+fs2))*fd/180*pi 
  1422.             aline f1*cos(f2) -f1*sin(f2) 
  1423.             f1 = fr-PhotonA*fs9 
  1424.             f2 = (1.45+GluonW*(1+fs1))*fd/180*pi 
  1425.             aline f1*cos(f2) -f1*sin(f2) 
  1426.             f1 = fr-PhotonA 
  1427.             f2 = (1.5+GluonW)*fd/180*pi 
  1428.             aline f1*cos(f2) -f1*sin(f2) 
  1429.             f1 = fr-PhotonA*fs9 
  1430.             f2 = (1.55+GluonW*(1-fs1))*fd/180*pi 
  1431.             aline f1*cos(f2) -f1*sin(f2) 
  1432.             f1 = fr-PhotonA*fs8 
  1433.             f2 = (1.6+GluonW*(1-fs2))*fd/180*pi 
  1434.             aline f1*cos(f2) -f1*sin(f2) 
  1435.             f1 = fr-PhotonA*fs7 
  1436.             f2 = (1.65+GluonW*(1-fs3))*fd/180*pi 
  1437.             aline f1*cos(f2) -f1*sin(f2) 
  1438.             f1 = fr-PhotonA*fs6 
  1439.             f2 = (1.7+GluonW*(1-fs4))*fd/180*pi 
  1440.             aline f1*cos(f2) -f1*sin(f2) 
  1441.             f1 = fr-PhotonA*fs5 
  1442.             f2 = (1.75+GluonW*(1-fs5))*fd/180*pi 
  1443.             aline f1*cos(f2) -f1*sin(f2) 
  1444.             f1 = fr-PhotonA*fs4 
  1445.             f2 = (1.8+GluonW*(1-fs6))*fd/180*pi 
  1446.             aline f1*cos(f2) -f1*sin(f2) 
  1447.             f1 = fr-PhotonA*fs3 
  1448.             f2 = (1.85+GluonW*(1-fs7))*fd/180*pi 
  1449.             aline f1*cos(f2) -f1*sin(f2) 
  1450.             f1 = fr-PhotonA*fs2 
  1451.             f2 = (1.9+GluonW*(1-fs8))*fd/180*pi 
  1452.             aline f1*cos(f2) -f1*sin(f2) 
  1453.             f1 = fr-PhotonA*fs1 
  1454.             f2 = (1.95+GluonW*(1-fs9))*fd/180*pi 
  1455.             aline f1*cos(f2) -f1*sin(f2) 
  1456.             f1 = fr 
  1457.             f2 = 2*fd/180*pi 
  1458.             aline f1*cos(f2) -f1*sin(f2) 
  1459.             end rotate 
  1460.         next fi 
  1461.         begin rotate fa-2*fn*fd 
  1462.             amove fr 0 
  1463.             f1 = fr+PhotonA*fs1 
  1464.             f2 = (0.05+GluonW*(1-fs9))*fd/180*pi 
  1465.             aline f1*cos(f2) -f1*sin(f2) 
  1466.             f1 = fr+PhotonA*fs2 
  1467.             f2 = (0.1+GluonW*(1-fs8))*fd/180*pi 
  1468.             aline f1*cos(f2) -f1*sin(f2) 
  1469.             f1 = fr+PhotonA*fs3 
  1470.             f2 = (0.15+GluonW*(1-fs7))*fd/180*pi 
  1471.             aline f1*cos(f2) -f1*sin(f2) 
  1472.             f1 = fr+PhotonA*fs4 
  1473.             f2 = (0.2+GluonW*(1-fs6))*fd/180*pi 
  1474.             aline f1*cos(f2) -f1*sin(f2) 
  1475.             f1 = fr+PhotonA*fs5 
  1476.             f2 = (0.25+GluonW*(1-fs5))*fd/180*pi 
  1477.             aline f1*cos(f2) -f1*sin(f2) 
  1478.             f1 = fr+PhotonA*fs6 
  1479.             f2 = (0.3+GluonW*(1-fs4))*fd/180*pi 
  1480.             aline f1*cos(f2) -f1*sin(f2) 
  1481.             f1 = fr+PhotonA*fs7 
  1482.             f2 = (0.35+GluonW*(1-fs3))*fd/180*pi 
  1483.             aline f1*cos(f2) -f1*sin(f2) 
  1484.             f1 = fr+PhotonA*fs8 
  1485.             f2 = (0.4+GluonW*(1-fs2))*fd/180*pi 
  1486.             aline f1*cos(f2) -f1*sin(f2) 
  1487.             f1 = fr+PhotonA*fs9 
  1488.             f2 = (0.45+GluonW*(1-fs1))*fd/180*pi 
  1489.             aline f1*cos(f2) -f1*sin(f2) 
  1490.             f1 = fr+PhotonA 
  1491.             f2 = (0.5+GluonW)*fd/180*pi 
  1492.             aline f1*cos(f2) -f1*sin(f2) 
  1493.             f1 = fr+PhotonA*fs9 
  1494.             f2 = (0.55+GluonW*(1+fs1))*fd/180*pi 
  1495.             aline f1*cos(f2) -f1*sin(f2) 
  1496.             f1 = fr+PhotonA*fs8 
  1497.             f2 = (0.6+GluonW*(1+fs2))*fd/180*pi 
  1498.             aline f1*cos(f2) -f1*sin(f2) 
  1499.             f1 = fr+PhotonA*fs7 
  1500.             f2 = (0.65+GluonW*(1+fs3))*fd/180*pi 
  1501.             aline f1*cos(f2) -f1*sin(f2) 
  1502.             f1 = fr+PhotonA*fs6 
  1503.             f2 = (0.7+GluonW*(1+fs4))*fd/180*pi 
  1504.             aline f1*cos(f2) -f1*sin(f2) 
  1505.             f1 = fr+PhotonA*fs5 
  1506.             f2 = (0.75+GluonW*(1+fs5))*fd/180*pi 
  1507.             aline f1*cos(f2) -f1*sin(f2) 
  1508.             f1 = fr+PhotonA*fs4 
  1509.             f2 = (0.8+GluonW*(1+fs6))*fd/180*pi 
  1510.             aline f1*cos(f2) -f1*sin(f2) 
  1511.             f1 = fr+PhotonA*fs3 
  1512.             f2 = (0.85+GluonW*(1+fs7))*fd/180*pi 
  1513.             aline f1*cos(f2) -f1*sin(f2) 
  1514.             f1 = fr+PhotonA*fs2 
  1515.             f2 = (0.9+GluonW*(1+fs8))*fd/180*pi 
  1516.             aline f1*cos(f2) -f1*sin(f2) 
  1517.             f1 = fr+PhotonA*fs1 
  1518.             f2 = (0.95+GluonW*(1+fs9))*fd/180*pi 
  1519.             aline f1*cos(f2) -f1*sin(f2) 
  1520.             f1 = fr 
  1521.             f2 = (1+2*GluonW)*fd/180*pi 
  1522.             aline f1*cos(f2) -f1*sin(f2) 
  1523.         end rotate 
  1524.             end if 
  1525.         end origin 
  1526.     amove ffx ffy 
  1527.     end if 
  1528. end sub 
  1529.  
  1530. ArrowL = 0.5  ! length 
  1531. ArrowA = 0.2  ! amplitude 
  1532. ArrowD = 0    ! displacement 
  1533.  
  1534. ! Arrow on the last line, fraction F from the beginning 
  1535. sub Arrow f 
  1536.     if fl>0 then 
  1537.         gsave 
  1538.         amove fx fy 
  1539.         begin origin 
  1540.             if fr>0 then 
  1541.                 f1 = fa+fc*(f*fl+90) 
  1542.                 fx1 = 0 
  1543.                 fy1 = ArrowD-fc*fr 
  1544.             else 
  1545.                 f1 = fa 
  1546.                 fx1 = f*fl 
  1547.                 fy1 = ArrowD 
  1548.             end if 
  1549.             begin rotate f1 
  1550.                 amove fx1 fy1 
  1551.                 begin origin 
  1552.                     aline -ArrowL ArrowA 
  1553.                     amove 0 0 
  1554.                     aline -ArrowL -ArrowA 
  1555.                 end origin 
  1556.             end rotate 
  1557.         end origin 
  1558.         grestore 
  1559.     end if 
  1560. end sub 
  1561.  
  1562. MomL = 0.75 ! half-length 
  1563. MomD = 0.5  ! displacement 
  1564.  
  1565. ! Arrow for indication of momentum near the last line 
  1566. sub Mom f 
  1567.     if fl>0 then 
  1568.         gsave 
  1569.         amove fx fy 
  1570.         begin origin 
  1571.             if fr>0 then 
  1572.                 f1 = fa+fc*(f*fl+90) 
  1573.                 fx1 = 0 
  1574.                 fy1 = MomD-fc*fr 
  1575.             else 
  1576.                 f1 = fa 
  1577.                 fx1 = f*fl 
  1578.                 fy1 = MomD 
  1579.             end if 
  1580.             begin rotate f1 
  1581.                 amove fx1 fy1 
  1582.                 begin origin 
  1583.                     amove MomL 0 
  1584.                     aline -MomL 0 
  1585.                     amove MomL 0 
  1586.                     aline MomL-ArrowL ArrowA 
  1587.                     amove MomL 0 
  1588.                     aline MomL-ArrowL -ArrowA 
  1589.                 end origin 
  1590.             end rotate 
  1591.         end origin 
  1592.         grestore 
  1593.     end if 
  1594. end sub 
  1595.  
  1596. ! Double arrow near the last line 
  1597. sub DMom f 
  1598.     if fl>0 then 
  1599.         gsave 
  1600.         amove fx fy 
  1601.         begin origin 
  1602.             if fr>0 then 
  1603.                 f1 = fa+fc*(f*fl+90) 
  1604.                 fx1 = 0 
  1605.                 fy1 = MomD-fc*fr 
  1606.             else 
  1607.                 f1 = fa 
  1608.                 fx1 = f*fl 
  1609.                 fy1 = MomD 
  1610.             end if 
  1611.             begin rotate f1 
  1612.                 amove fx1 fy1 
  1613.                 begin origin 
  1614.                     amove MomL*(1-(DoubleA/ArrowA)) DoubleA 
  1615.                     aline -MomL DoubleA 
  1616.                     amove MomL*(1-(DoubleA/ArrowA)) -DoubleA 
  1617.                     aline -MomL -DoubleA 
  1618.                     amove MomL 0 
  1619.                     aline MomL-ArrowL ArrowA 
  1620.                     amove MomL 0 
  1621.                     aline MomL-ArrowL -ArrowA 
  1622.                 end origin 
  1623.             end rotate 
  1624.         end origin 
  1625.         grestore 
  1626.     end if 
  1627. end sub